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Abstract 

I present an algorithm based on sector decomposition and Mellin-Barnes techniques to power 
expand Feynman integrals. The coefficients of this expansion are given in terms of finite integrals 
that can be calculated numerically. I show in an example the benefit of this method for getting 
the full analytic power expansion from differential equations by providing the correct ansatz for 
the solution. For method of regions the presented algorithm provides a numerical check, which is 
independent from any power counting argument. 
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I. INTRODUCTION 



For power expanding Feynman integrals several methods exist, where all of them have 
their limitations. Mellin-Barnes techniques provides a very general method to obtain all 
powers 111 y|- This method however fails if the integrals are getting too complex. On the 

nnnn 

other hand method of regions [2|, y, U, [5( is a convenient way to obtain the leading power, 
whereas it is getting rather complicated for higher powers because of the many contributing 
regions and because it is difficult to automatize. Furthermore it is a very non-trivial task to 
make sure that one has not forgotten or counted twice any region. However in the Euclidean 
imit, where no collinear divergences arise, automatizations exist, which rely on graph theory 
f], Ql. Another way to expand Feynman integrals, which has been proposed and worked out 



in 



111 ], is based on differential equations. Differential equation techniques, which 



has been proposed first in 12j, is easy to automatize in a computer algebra system. This 
makes it a convenient method to obtain subleading powers, whereas the leading power is 
in most cases needed as an input like a boundary condition. Another limitation is the fact 
that this method relies on a correct ansatz in terms of powers of the expansion parameter. 
However it is a priori not obvious which powers of the expansion parameter occur (e.g. only 
integer powers or also half-integer powers). 

In the present paper I present a semi-numerical method, that provides the power ex- 
pansion of Feynman integrals by giving explicit expressions of the expansion coefficients in 
form of finite integrals the can be solved numerically. In particular this method gives the 
contributing powers of the expansion parameter, from where one can read off the correct 
ansatz to solve the differential equations that determine the set of Feynman integrals. 

The algorithm that is worked out in the present paper combines sector decomposition 
[3, 14, [3 [3] with Mellin-Barnes techniques. It is completely independent from any power 
counting argument such that it can be used as a cross check for method of regions. This is 
very useful in cases, where method of regions becomes involved because of many contributing 
regions. 

The paper is organized as follows. In Section [IT] the algorithm is explained in detail. 
In Section II III I apply this algorithm to a set of two Feynman integrals, that are power 
expanded by differential equation techniques, where the leading powers are obtained by 
method of regions. I will show explicitly how this algorithms gives the correct ansatz for 
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the differential equations and provides a non-trivial check for method of regions. 



II. ALGORITHM 



We follow the steps of Section 2 of [13j|. We start with a L-Loop Feynman integral 



l(2n) D P 1 ...P N 
which using the Feynman parameterization 

*(i-Eii* 



(i) 



- 1 - =T(N) f d N x 
Pi---Pn Jo (x 1 P 1 + ... + x n P n ) n 



(2) 



can be cast into the form: 



- N L 

G = T(N) / d N x5(l-J2^n) / I] 

n=l J i=l 



d D ki 

D 



A (27T) 



^2 k j ■ hMji - 2 V" kj ■ Qj + J 

3,1=1 3=1 



-N 



(3) 



We define D = 4 — 2e as usual. After performing the integration over the loop momenta we 
obtain: 



G=(-l) 



N 



(4tt)^/ 2 



N 



T(N — LD/2) / d N x5(l-J2 



n=l 



jjN-(L+1)D/2 
' pN-LD/2 ' 



where 



and 



F = -det(M) 



3,1=1 



(4) 



(5) 



U = det(M). (6) 

Let us assume © contains the parameter A, in which we want to expand ©• Using the 
Mellin-Barnes representation 

1 1 1 rio ° 



/lOO 
dsT(-s)T(s + x)X°X, 
■ioo 



(Xi + X 2 y T(x) 2m 
where the integration contour over s has to be chosen such that 



2 ) 



(7) 



-x < K(s) < 0, 
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we modify in the following way 

G=(-V N (ttTTTT? ) 7T-. I ds\ s T(-s)T(s + N-LD/2) 



\ L 1 
i \ 1 



f OC 



x / d N x 5(1 - ^ x n )U N -^ D / 2 F°F 2 s - N+LD/ \ (8) 

^ n=l 



where 



F = XFi + F 2 . (9) 

The main idea behind the procedure below is the following: By closing the integration path 
to the right hand side of the imaginary axis we sum up all the residua on the positive real 
axis and obtain an expansion in A. Powers of In A appear because of poles of order higher 
than one and because of terms of the form X A ~ Bt in the expansion in A. These terms turn 
after expanding in e into powers of In A. 

We continue with part I and II of [3]. First we split the integral over the Feynman 
parameters into 

» n „ N 

d N x = J2 d N x II 6 ( x * ~ x ^ ( 10 ) 

J 1=1 3=1 

and integrate out the 5-function by the substitution 

X t tj j < I 

xi j = l (11) 
Xitj-i j > I 



such that we obtain 

G = (-1) N (jjjyw) L ^ [2 ds m-s)r(s + N- LD/2) £ jf * d»~H G h (12) 
where 

G l = U»-W°/* F S lF - r N+LD/2 (13) 

is obtained by the substitution (iTTj) . In (|T2|) the integration over small t leads to poles in s. 
This behavior is made explicit, if we follow the steps of Part II of jl^]: Look for a minimal 
set {t ai , . . . ,t ar } such that Ui, Fij or F 2> i vanish, if these parameters are set to zero. We 
decompose the integral into r subsectors 



fd N -H= /v^n^-t. 
j ° j ° k =i i= i 



(14) 



and substitute 

t aj ^\ takta i ^\ (15) 

which leads to the Jacobian factor t r ~ k l . Now we are able to factorize out t ak from Ui, F\j 
or F 2 j. After repeating these steps, until Ui, F^i and F 2 j contain terms that are constant 
in t, we end up with integrals over the Feynman parameters of the form 

EE / dN ' H II Ut {L+ ^ ,2 FU k F-- N ^\ (16) 

1=1 k Jo \j=i J 

where Ufa, Fi t fa and F 2y ik contain terms that are constant in t. T he p rocedure above can in 
principal lead to infinite loops. This problem was addressed in lj[l^], where algorithms 
are proposed that avoid these endless loops by choosing appropriate subsectors. I have not 
yet faced any endless loop in the problems I dealt with. However one should keep in mind 
that they can occur and adapt the implementation of the algorithm if needed. 
From ( jTBT) we can read off that the poles in s are located at: 

s^ 1 + n+ £- B > C , (IT) 

where n G No- Eq. ( |T7|) becomes clear if one Taylor expands in ( |T6l) the terms outside the 
brackets with respect to tj and performs the integration. 

In f ll2p we have to choose the contour of the integration over s such that the integration 
over the Feynman parameters tj converges. This leads to the condition 

Aj - Bje - Cj^(s) > -1 Vj. (18) 

The poles in ( 1171) that have to be taken into account are those that are located on the right 
hand side of the integration contour, i.e. 

< s jn . (19) 

From (JIZD and ([IE} we conclude that (EE]) is fulfilled if and only if Cj > 0. 

In the next step we calculate the residue of ( jTBT) at sj n . We write the fc'th Feynman 
integral in the form 

1 dt k t*~ B ' e ~ C ' {S ~ Si Mh,s) (20) 

o 
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and note that this term is singular in s — Sj n if and only if 

B' = and A' < -1. (21) 

n 

So following Part III of [13j we expand I(tk, s) around t k = and obtain 

-A'-l P 

T{t k ,s)= Y,lto\s)± + R(t k ,s), (22) 
with a rest term R(tk, s) = 0(t~ A '), such that ( 1201) becomes 

where we used that B' = 0. We repeat this procedure for all fc where condition ( |2T|) is 
fulfilled. The remaining integrals do not diverge for s = Sj n . So it is save to expand them 
around s — Sj n and we can easily calculate the residue at s — Sj n . 

What is left is to calculate the Laurent expansion in e. From the previous procedure we 
obtain terms of the form 

J\ n t (j{tf' B ^{\nt^) I(t,e) (24) 

- C '{s-s jn ) 

The logarithms (lnt,-)^ arise from taking the residues of terms of the form ^ m with 

m > 2. In (124j) we wrote these logarithms explicitly such that we can expand I(t, e) around 
tj = 0. The poles in e in (124p originate from integrals 

dtjtf-^QntjpIfae) (25) 
with A" < —1. Repeating the procedure above we expand 

— A" — 1 

X{t jt e)= ^\ed +R{t jl e) (26) 

and obtain for ( 1251 ) 



' A'l-B'Je. 



^ (1+P + A'j - B'Je) a i +1 p\ Jo 3 

All the remaining integrals over tj are finite and can in principle be calculated numerically. 
Finally the original integral G in ([3]) obtains the form 

G= ^e l A m (lnA)% m , n , (28) 

i,m,n 



k + l 



FIG. 1: Sunrise diagrams. The thick line denotes a propagator of mass M, while the thin lines 
stand for mass m. The double line denotes that the propagator is to be taken squared. 

where the Ij iTnjTV contain finite integrals that can be numerically evaluated. The logarithms 
(In A) n arise both due to poles of higher order in the Mellin-Barnes parameter and to the 
expansion in e from terms of the form X € /e n . Depending on the values of Cj in (I17p the sum 
over m does not only run over integer numbers but also over numbers of the form 

1 + n + Aj 

where n is integer. I stress that even if a numerical evaluation of the integrals i» mn is not 
possible, we can obtain non-trivial statements about the power expansion of G from (1281) 
together with (TP?]) . That is to say ( fTTl) gives us information about the possible powers of A 
e.g. we know if we only get integer powers or also powers of And from (128]) we can read 
off up to which power In A appears. As we will see in the next section this information will 
prove to be useful to obtain the power expansion by means of differential equations. 



III. EXAMPLE: POWER EXPANSION OF FEYNMAN INTEGRALS BY DIF- 
FERENTIAL EQUATION TECHNIQUES 

The idea to get the expansion of Feynman integrals by differential equations has been 
proposed and worked out in , lUj. By the following example we will see that the 



algorithm shown in the last section will give us the correct ansatz to solve the given system 
of differential equations and help us with the calculation of the initial conditions. We start 
with the integrals given by Fig. [T], where we assume p 2 = M 2 : 

d D k d D l 1 



h 



(2vr) D (2n) D (k 2 + 2k ■ p) ((k + I) 2 - m 2 ) (Z 2 - m 2 ) 
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d^J^l 1 m 

(2n) D (2n) D (k 2 + 2k-p){{k + l) 2 -m 2 ){l 2 -m 2 ) 2 ' V ; 

Let us assume that we want to expand these integrals in A = m 2 /M 2 and need the result 

up to order e. For simplicity let us also set M 2 = 1 and m 2 = X. Using integration-by-parts 



identities 



19 



20 



with 



and 



h 



2l(], we get the following differential equations for 1\ and l 2 . 

^-h = huh + h 12 I 2 + 9i 
dX 

d 

—I 2 = h 2 ih + h 22 I 2 + 92 
dX 



2 

1-3A 

2A(1-A) 2A(1-A) 




1 



+ e 





-2+4A 

4A(1-A) A(l-A) 
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2A(1-A) 



(30) 



(31) 



92 



'l-e] 



d D k d D l 



1 



4A 2 (1 - A) [J (2n) D (2n) D ((Jfe + I) 2 - A) (I 2 - A) 

d D k d D l 1 " 

(2tt) d (2ir) D (k 2 + 2k-p) (I 2 - A) 

i \l-2e \-e 

: rr e ) 2A ~ A 



(4rr) D v ^ 4A(1-A) 
where ( l3~TT) and (1321) are exact in A and e. By defining 

= E^f^QnA)' 



9* = E^f^i^ 



301) becomes 

o = o-+i)i;;:r M + a- + i)/. 



2 j 

(j'+i,*) i ^ i ^ r 0'+l,W)_ Y" rO'-j".*) _ _0'.fc) 

a,i (l ai3,i' 1 p,i-i' 9a,i ■ 

f3=l,2i'=0 j'=-l 



(32) 



(33) 



(34) 



In ( 1331 we have not yet specified which values the summation index j takes and up to which 
maximum value the finite sum over k runs. By implementing the steps of the last section, 
which led to (1171) . in a computer algebra system we obtain from (1T71) that Ji comes with the 
powers of A 

(35) 



A n , X n+1 -\ 



n + 3 



A 2 



-2c 



and I 2 with 



A n , \ n ~\ A^~ 2e , (36) 



where n G No- From ( !35l) and ( |36l) we read off that j takes the values 0, 1/2, 1, . . .. In (JB! 
integer-valued and half-integer-valued j do not mix. So we would have missed powers of 
vA, if we had made the naive ansatz that only come with integer powers of A. Now 
one could argue that vA is already contained in the sum over In A. However in order to 
solve we have to assume that there exists /c max such that I^f' = for all k > /c max - 
A computer algebra analysis of the algorithm in the previous section tells us that in our 
special case fc max = 3. 

Solving ( 1341) up to (9(e) we note that we need l[°f^ and /j 2 / ' as initial conditions, which 



can be obtained by method of regions 



3, Q, E 



5]. In the case of we note that only the 



region participates where both integration momenta are hard: 

k" = 0{l) and P = C(1). (37) 



In this region we obtain 

d D k d D l 1 1 T(-l + 2e)r(e)r(l -e) 2 r(3-4e) 

{2tt) d {2tt) d {k 2 + 2k-p){k + l) 2 l 2 ~ {4n) D r(2 - 2e)r(3 - 3e) : 



(38) 



(1,0) 

which is the leading power of I±. For I 2 Y we need the region where both k and / are soft, 
i.e. 

p = C( v / A) and F = 0(V\). (39) 
This region starts participating at 0{\f\): 

d D k d D l 1 

(2ir) D (2n) D (2k ■ p) ((k + /) 2 - A) (I 2 - A) 2 ~ 

-1 2- 2e 7rr (e - i) T (2e - |) i „ 

V 2 ) V ?J-\2- 2 \ (40) 



(An) D r(e) 

By comparing these results to flHol) and flHBT) we note that flHHl) and fHUj) correspond to definite 
poles in the Mellin-Barnes representation i.e. at s = and s = 1/2 — 2e. By (1171) and ( |2"3"|) 
we can calculate the coefficients of A and \^~ 2e in the A-expansion of I± and -Z2 numerically. 
This is a non-trivial test that we have not forgotten a contributing region, which is in general 
a problem of method of regions. 
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We normalize our integrals by multiplication with (exp(7E)/(47r)) 2,E and obtain from the 
solution of (Ellj) the analytical expansion in A and e: 



h 



(47T) 



A 



2e 2 4e 
1 -3 + 21nA 

5 2 e 
3vr 2 26 



5 11 5tt 2 /55 25tt 2 11 , 

el-- — - ^rC(3) + 



12 



7T 



16 24 



+ y + 6 In A- (lnA) 2 + 



7T 



M "3 + ~— + fC(3) + I 14 + - \ In A - 3(lnA) 2 + , 



(lnA) ; 



, 3 — 167T 

A^e — - — + C(A 2 ) 



+ 0(e 2 ) 



+ 



(4tt) 



1 -1 + 2 In A 1 7T 2 , , 1, , x2 
— ? + + - + — + 2 In A - -(In A) + 

2e 2 2e 2 4 2 K J 

11 U7r2 + yC(3) + (4 + In A - (In A) 2 + i(lnA) 3 ) + 



2 12 
A^(-4e7r 2 ) + 



A 



l-y + lnA-^(lnA) 2 + 
e ( 11 + ^ - 4C(3) - 3 In A - i(ln A) 2 + ~(ln A) 3 



+ 



A f e lj! + 0(A 2 ) 



G(e 2 



On the other hand our numeric method of Section [TT] gives 
h = 





0.20 


10- 4 




e 2 


*(- 


0.40 


— + 



+ 



-2.2-4.5e+ 



- 0.026 + 2.4 In A - 0.40(ln A) 2 + 



e (8.9 + 6.9 In A - 1.2(lnA) 2 + 0.131nA) 3 ) J - 



21.eA2 + 0{\ 2 



+ 0(e 2 



(41) 



10 



10 



,-4 



0.20 



0.20 + 0.40 In A 



+ 1.2 + 0.80 In A - 0.20(lnA) 2 + 



+ 



e 



e (7.9 + 2.2 In A - 0.40(ln A) 2 + 0.067(ln A) 3 ) 



16.eA2 + 



A y - 1.7 + 0.40 In A - 0.20(ln A) 2 + 

e(5.1 - 1.21nA-0.20(lnA) 2 + 0.20(hiA) 3 )^ + 
5.3eAi + C(A 2 ) + £>(e 2 ), 



(42) 



which is consistent with (l4"Ii) . 

IV. CONCLUSIONS 

By combining sector decomposition with Mellin-Barnes techniques I developed an algo- 
rithm for power expanding Feynman integrals, where the coefficients in the expansion are 
given by finite integrals. Even if these integrals cannot be evaluated numerically, we can 
read off, which powers of the expansion parameter contribute and up to which power the 
logarithms occur. This non-trivial information provides the correct ansatz for solving the 
set of differential equations that determine the Feynman integrals. 

Another application of the presented algorithm is testing method of regions numerically. 
We have seen that every region, that has a unique scaling in the expansion parameter, cor- 
responds to a definite power in the Mellin-Barnes expansion. So it can be tested separately. 
For method of regions it is often an involved problem to make sure not to have missed 
or counted twice any region. This algorithm provides a test of method of regions that is 
independent of any power counting argument. 
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